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Abstract 


In this report analytical and numerical methods are proposed to estimate the effective 
elastic properties of regular and random open-cell foams. The methods are based on the 
principle of minimum energy and on structural beam models. The analytical solutions are 
obtained using symbolic processing software. The microstructure of the random foam is 
simulated using Voronoi tessellation together with a rate-dependent random close- 
packing algorithm. The statistics of the geometrical properties of random foams 
corresponding to different packing fractions have been studied. The effects of the packing 
fraction on elastic properties of the foams have been investigated by decomposing the 
compliance into bending and axial compliance components. It is shown that the bending 
compliance increases and the axial compliance decreases when the packing fraction 
increases. 
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X. Introduction 


Structural foams are made with high stiffness-to-weight ratio and high-energy 
absorption and are very suitable for lightweight structures and insulating materials [1-4]. 
Different materials can be foamed by a variety of manufacturing processes. For example, 
metals can be foamed by injecting gases or by adding gas-releasing blowing agents [1]. 
The foams produced exhibit different microstructures and may have open and closed cells 
with different degrees of randomness. 

Foams are highly heterogeneous materials with mechanical and acoustic 
properties determined not only by the frame material but also to a significant degree by 
their microstructure. To model the properties of these materials, one first needs to 
generate the microstructures. The actual foam microstructure can be obtained by 3-D X- 
ray tomography [5]. The microstructures can also be simulated. In analysis, foams are 
often approximately represented as a periodic array of regular cells such as hexagonal 
and truncated octahedron. The random foams may also be simulated with Voronoi 
tessellation [6-10]. In this approach the randomness of the foams is determined by the 
distribution of the “seed points” which are usually generated using hard spherical (disks) 
packing algorithms [10]. Another approach which is particularly useful in modeling 
bicontinuous media is the level-cut Gaussian random field scheme [10,1 1]. 

The elastic properties of the foams are closely related to their microstructures. For 
regular cells, by analyzing the smallest unit cell and representing each strut by a 
structural beam, analytical formulations were obtained for some 2-D and 3-D structures 
[12-17]. In these studies, the regular foams were usually assumed to be formed by cells 
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arranged in periodic arrays. Using thin beam or plate theories, the elastic constants 
(modulus and Poisson ratio) can be explicitly determined as functions of the 
microstructure of the cell. Using this approach, an analytic formulation for open cells 
with tetrakaidecahedral microstructure has been obtained [16, 17]. It has also been shown 
that the elastic properties for 3-D regular structures can be obtained by averaging in 
different directions [15] the properties obtained for 2-D structures. Empirical data shows 
that, for random open cell structures, the modulus is a power function of the relative 
density [2]. Rigorous analyses of random structures have been obtained using finite 
element analysis based on structural beam elements [7] and solid cubic elements [9]. 

In this report we describe a general approach to model the elastic properties of 
regular open-cell foams numerically and analytically, based on the principle of minimum 
energy and on structural beam models. Examples for 2-D general honeycomb, triangular 
and rectangular cells and 3-D orthogonal, rhombic dodecahedron and truncated 
octahedron are given. We found the cell compliance matrix can usually be decomposed 
into two terms which correspond to bending and axial compliances respectively. Also 
foams with different degrees of randomness have been simulated using a rate-dependent 
random packing algorithm and Voronoi tessellation. Then the effect of packing fraction 
on the elastic properties has been addressed. Here we will discuss the evolution of the 
geometry property of the random foams generated from the Voronoi diagram for both 2D 
and 3D using the rate-dependent algorithm described above. The corresponding elastic 
properties will be discussed in the following sections. 
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2. Simulation of 2D and 3D foam microstructures 


In order to understand the relation between foam microstructure (distribution of 
cell shape, size and flaws in the cells) and elastic properties, one needs to simulate open 
foam microstructures with specific geometrical properties. Due to the similarity between 
the Voronoi diagram and foam microstructures, the Voronoi tessellation has been widely 
used to simulate the foam structures [6-10], Here we will describe a procedure to 
generate foams with different regularities by combining a rate-dependent closed random 
packing algorithm with Voronoi tessellation. 

2.1. Voronoi tessellations 

In this simulation a random set of seed points or nuclei is first generated. Then 
the Voronoi tessellation algorithm is used to partition the space into cells, each of which 
consists of the area enclosing one particular nucleus (we have used the “Quickhull” [18] 
algorithm to produce Voronoi tessellation for 2D and 3D dimensions). For open-cell 
foams, the struts are formed by the edges of the Voronoi cells. 

The regularity of the generated cells is determined by the distribution of the 
nuclei. For regular distributions, foams with ordered uniform cells are generated. Figures 
1 and 2 show simulated regular 2-D and 3-D Voronoi cells whose elastic properties will 
be discussed in the following sections. Figure i shows the (a) rectangular, (b) triangular 
and (c) general hexagonal 2-D structures corresponding to rectangular, hexagonal and 
triangular nuclei respectively. The smallest repeating cells are shown in the bottom of the 
figures. Figure 2 shows 3-D structures: (a) cubic, (b) rhombic dodecahedron and (c) 
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truncated octahedron structures. The Voronoi unit cells are shown in the bottom of these 


figures. 


2.2. Random Voronoi tessellations 

The advantage of the Voronoi method is that it allows simulating microstructures 
with different level of randomness, which reflects the natural variations in microstructure 
in practical foams. The regularity of the simulated foams is determined by the 
randomness of the seed nuclei. There are different point-process algorithms to generate 
spatial points with different degrees of randomness [19]. The simple sequential inhibition 
and the close random packing point process can be used to produce a set of points widely 
varying in regularity. The Gibbs point processes are good models for structures with 
some degree of regularity and the Poisson point processes for random structures without 
regularity. 

Here we use the close random packing point process to generate spatial points 
(seed nuclei) with different randomness. The randomness of the spatial points is 
represented by the packing fraction (j>. Different computer algorithms have been proposed 
to achieve the maximum packing fraction <j> c . The maximum packing fraction <j> c of a rate- 
dependent densification algorithm [20, 21] can reach 0.649, of a Monte Carlo scheme 
[22] can achieve 0.68 and a “drop and roll” algorithm [23] 0.60. Because of the efficiency 
of the rate-dependent densification algorithm, we will apply it to generate random nuclei 
with different regularities and combine it with Voronoi tessellation to obtain foam 
structures with different regularities. In this section statistics of geometrical 
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characteristics of the cell are studied. The statistical description of the foam elastic 
properties will be discussed in the following sections. 

In the rate-dependent densification algorithm [20, 21], first a set of N full random 
points are generated in a rectangle (2D) or cubic (3D) box with edge length B with 
periodic boundaries as the starting state for iteration. Each point is the center of an inner 
and an outer circle for 2D or sphere for 3D. At iteration step /, the diameter di of the inner 
sphere or circle is the minimum distance between any two points. The corresponding 

minimum packing fraction is ^ for 3D or <j>, = f° r 2D. In each 


step, the two points with minimum distance are spread apart symmetrically along the line 
joining their centers until the distance between these centers are equal the outer diameter 
Dj. Then the diameter D, of the outer sphere or circle is reduced based on the minimum 


TtW ( D V 

packing fraction <j>; and the maximum packing fraction O, (<J>, = -^—1 I for 3D 


or 


nNf D V 

= TT for 2D) at step f [21i 


d m = A -£(*, -♦,)*, 


(i) 


where R and a are parameters to control the convergence rate. The initial value of the 

outer diameter is D 0 = B\ j for 3D and Z) 0 -B \— -| for 2D, which results in 

\nN ) l tc N ) 


maximum packing fraction <J> 0 = 1 . Iteration continues by spreading the new pair of 
points with minimum distance until the inner diameter is no less than the outer diameter. 
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As an example. Figure 3 shows the evolution of the minimum packing fraction. 
For random foam simulation, in the2D case we use 1000 nuclei points in a unit rectangle 
area with control parameters /?=0.00015 and a=0.01 and reach the final packing fraction 
0.816. For 3D foams, we use 2000 nuclei points in a unit cubic volume with control 
parameters i?=0.0001 and a=0.3 and reach the final packing fraction is 0.643. Although 
the packing fraction is not monotonically increasing with iteration number, in general the 
packing fraction steadily increases and the variation amplitude decreases. The point sets 
corresponding to the iteration on the upper profile represent points with different packing 
fraction. Using the point sets obtained with different packing fraction, the corresponding 
Voronoi diagrams are obtained by the “Quickhull” algorithm. Examples for 2D and 3D 
foam structure evolutions with increase of packing fraction are shown in Figures 4 and 5. 
Zero packing fraction corresponds to the initial fully random state. As can be seen, as 
packing fraction increases, the size and shape of the cells become more uniform. 
Statistics of the geometrical properties of the random foams will be discussed in the next 
section. 

2.3. Distribution of geometrical parameters of random foams 

Zhu et al. [24] use a “drop and roll” algorithm [23] to generate random nuclei and 
studied the geometric properties of the corresponding Voronoi diagram. However the 
“drop” algorithm can only achieve packing fracture 0.5 for 2D foams. Jullien et al [21] 
used the rate-dependent algorithm to study the long-range correlation and local order in 
random packing of spheres via calculation of the two-point correlation function and the 
Voronoi tessellation. Using a discrete element method, Yang et al. [25] studied the 
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topological and metric properties of the 3D Voronoi diagram as a function of the packing 
particle size and fraction. 

Here we will discuss the evolution of the geometric properties of the random 
foams generated from the Voronoi diagram for both 2D and 3D using the rate-dependent 
algorithm described above. Figure 6 reports the distribution functions of the cell 
geometrical parameters for a 2D random structure at different packing fractions. The 
geometrical properties (edge length, vertex angle, cell area and cell perimeter, number of 
edges) of the cells are normalized by the corresponding parameter of the regular 
honeycomb cell. The honeycomb cell corresponds to the highest packing fraction 



« 0.9069) for 2D structures. In Fig. 6(a), the distribution of edge length is 


shown. For fully random packing of nuclei, the normalized edge length is almost 
uniformly distributed between 0.0 and 1.5 and no obvious peak can be observed. As the 
packing fraction increases, a peak above 1 appears and approaches 1. Figure 6(b) shows 
the distribution of the cell perimeter. In contrast to the edge length shown in Fig. 6(a), the 
peak position for the cell perimeter is much closer to 1 and has no observable shift as the 
packing fraction increases. Figure 6(c) presents the distribution for vertex angle which is 
normalized by the vertex angle 120° for a regular honeycomb. For lower regularity, the 
distribution peak is above 1 and approaches 1 as regularity increases. Figure 6(d) shows 
the cell area distribution. As can be seen, the normalized cell area is less than 1 and 


approaches 1 as regularity increases. Figure 6(e) shows the distribution of edge numbers 
per cell. One can see that the peak is at 6 which is the edge number for a regular 
honeycomb. Only the probability at 6 increases as the foam regularity increases. All 
parameter distributions narrow with increase of the packing fraction and the peak position 
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is shifted to 1 if the cell parameter is normalized by the packing number of the 
honeycomb structure. These statistical distributions show that the full random structure 
evolves to regular honeycomb structures with number of iterations and increase of 
packing number 

Figure 7 shows the density ratio for the random foam generated from the Voronoi 
diagram. In this calculation, the Voronoi edges are replaced by struts. The width of the 
strut is selected to make the density ratio for the regular honeycomb structure 0.1. The 
relative density decreases linearly by 6% as packing fraction increases and can be written 
by fitting the curve to find the decrease rate 

— = 0.08767(<|> ( , -<|),) + 1 , (2) 

P* 


where p r and p A are the relative densities for the random and honeycomb structures 


respectively. <)>,,(= 0.9069) and <j) r are the packing fraction of the random and honeycomb 

2D structures respectively. The foam density decreases with packing. 

Figure 8 shows the distribution functions for 3D random open foam structures at 
different packing fractions. For the 3D case, the maximum packing fraction 


n 


(4> c = — 1 = = 0.74) is given by the face-centered cubic or hexagonal closed packing. 
3V2 


Their corresponding Voronoi diagram is the rhombic dodecahedron cell shown in Fig. 
2(b). However, as will be shown below, the final random close packing obtained using 
the procedure described above does not converge to the face-centered-cubic or the 
hexagonal-closed-packing structures. Foams with truncated octahedron cells have been 
widely studied [16, 17], therefore we will normalize the geometrical properties of the cell 
of the random structures by the parameters of the truncated octahedron, the body centered 
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packing structure (<j> c = = 0.68) whose cell is the truncated octahedron (Fig. 2(c)). 
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In Fig. 8(a), the distribution of edge length is shown. For full random packing of 
nuclei, the normalized edge length is dominated by short struts. With increase of packing 
regularity a peak at position larger than 1 starts to appear. As in the 2D case, very short 
edges still exist in highly regular structures. Figure 8(b) show the distribution of face 
areas. The area is normalized by the area of the hexagonal face in the truncated 
octahedron cell. Consistent with the edge length distribution, a significant number of 
faces have small area and the peak appears at slightly below 1. Figure 8(c) shows the 
distribution of vertex angle which is normalized by 120°, the vertex angle for the 
truncated octahedron cell. This distribution is very similar to the 2D case shown in Fig. 
6(b). The distributions of the surface area and volume of the Voronoi cells are shown in 
Figs. 8(c, d). Compared to Figs 8(a, b, c), these properties have much narrow 
distributions and sharper peaks. Figures 8 (f, g) show the distribution of number of facets 
per cell and number of edges per facet. The truncated octahedron has 14 facets per cell, 6 
of which are square and 8 hexagonal. The numbers corresponding to the truncated 
octahedron are indicated in Fig. 8(f, g) by the dotted lines. The peak for number of facets 
per cell is 15. The average number of facets per cell decreases from 15.35 to 14.25 when 
the packing fraction increases from 0 to 0.64. Figure 9 shows the density ratio versus 
packing fraction. In this calculation the Voronoi edges are replaced by thin struts. The 
width of the strut is selected to make the density ratio for the truncated octahedron foam 
10%. The density decreases by 10% as packing fraction increases. As in the 2D case, one 
may approximate the relation between the density ratio and packing fraction with a linear 
function: 
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— = 0.08767(ij), — <j) r ) -l- 1 , 

P, 


( 3 ) 


where p r and p, are the relative densities for the random and truncated octahedron 
structures respectively. 0.68) and i)> r are the packing fractions of the random and 
truncated octahedron 3D structures respectively. 

3. Analytical and Numerical Methods to Determine Elastic Properties of 
Foams 

Different approaches have been used for the homogenization of heterogeneous 
media and obtaining their effective elastic properties [10]. Here we will apply the 
principle of minimum potential energy [9, 10] to estimate the effective properties of 
foams. In this approach, we first simulate a foam structure bounded by a square for 2D 
and a cubic for 3D with periodic boundary conditions. Then all struts in the foam are 
modeled by plane or spatial beams using 2-nodes elements considering axial, bending 
and torsion deformations [26]. The strain energy of the system is calculated as the 
summation of the energy for each individual beam using the finite element method for 
spatial beams and the total strain energy U B is written as 

U b =~D t KD, (4) 

where D is the general displacement vector of the nodes which include axial, bending and 
torsion displacements; K is the global stiffness matrix. 

We consider that the macroscopic strain and stress o® are known at the 
computation domain boundary determining the boundary conditions. If assuming 
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macroscopic strain y®. (y® = e® for / = j and y® = 2s® for / * j ), the periodic boundary 
conditions are written as the relation between the displacements ( u i0 , u n , /'= 1, 2, 3) at the 
opposite side of the computation domain (square for 2D and cube for 3D) 

u io ~ u n =zl(x J0 -Xj 0 ), (5) 

where x i0 and x a are the coordinates of the boundary nodes at the opposite side. For 

given macroscopic strain, one may also use constrained boundary conditions. In this case, 
the displacement at the boundary is given as 

U i ~ X J > ( 6 ) 

where Xj is the coordinate of the boundary nodes. The bending moments at these nodes 
are assumed to be zero. For given macroscopic stresses a®, we calculate the force for the 
beam nodes lying on the boundary surfaces. For example, considering area A c with N 
nodes on the boundary surface, the force in each node is a°Ac/N and the moment is zero 

(due to constant stresses cr® ). Using different boundary conditions may have different 
results [27], 

After applying the boundary conditions the total strain energy (Eq. 4) of the 
system is obtained and is minimized using the fast conjugate gradient method. 
Alternatively one can also formulate the global finite element stiffness matrix for the 
system and solve the linear system of equations to obtain the equilibrium state whose 
corresponding strain energy represents the minimum values. 
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After the minimum strain energy is obtained, the effective elastic stiffnesses or 
compliances are calculated by comparing them with those of an effective homogeneous 
medium. The strain energy of a homogeneous medium for given strains s® is given by 


U = — p °C f°V 

U c. 2 'J ijkl^kl' t’ 


(7) 


where C iJU is the effective stiffness tensor; V e is the volume of the structure. For given 
constant stress a® , the strain energy of the effective homogeneous medium is given by 




( 8 ) 


where S iJU is the effective compliance tensor. Equating the strain energy obtained by the 


beam model ( U B =U e = U a ), we obtain the effective stiffness tensor for constant strain 

and the effective compliance tensor for constant stress. For example, for a 2-D problem 
the strain energy for the effective medium may be written as 

= K(( £ u) C m i + ( e 22) C 22 22 + 2E l ,e 22 C u2 2 +(Yi 2 ) C 12t2 )/2. (9) 

If Ejj is given and other strain components are zero, one obtains C 1U1 =2U t l(V e s° n 2 ) . 

Similarly one can obtain other constants. 

If the bounded structure includes a large number of struts, then the solution can 
only be obtained numerically. However for regular foams shown in Figs. 1 and 2, we 
have obtained simple analytical solution for the smallest cell in the foam using symbolic 
software. The results are discussed in the next two sections. 

For foams with multiple cell structures, one can first formulate the 
stiffness/compliance matrix for each cell and then obtain the elastic constants for the 
foams based on the connection of the cells. 
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4. Foams with 2-D Microstructures 

The elastic properties for 2D regular foam have been widely investigated [2, 15, 
14]. Here we use the method described in section 3 to drive the analytical elastic 
constants for 2-D foams. Differently from other approaches, we write the elastic 
constants in compliance form with shear deform correction. Therefore they could be used 
for forms with relatively larger thickness to length ratio or larger density ratio. We apply 
the algorithm described in section 3 for 2-D orthogonal, triangle and general honeycomb 
unit cells shown in Fig. 1. The analytic compliances for the unit cell shown in Fig. 1 are 
listed in Tables 1-3. All the compliances can usually be separated into two terms 
corresponding to bending and extension respectively. The entire cell compliance matrix 
Sfj can be written in the form as 


5 , a r 


( 10 ) 


where A is the cross-section area, I is the second moment of the cross section area and 
S s =1/ E S (E S is the modulus of the frame). Here we assume the second moments of area 

about the x and y axis are the same. The two parameters C and B depend on the geometry 
of the cell. The first term corresponds to axial compliance and the second to bending 
compliance. Tables 1-3 list the expressions for C and B for 2-D for orthogonal, triangle 
and general honeycomb cells. In these tables <t> y is the shear correction constant in the y 


axis [26]. The results obtained from the macroscopic strain or stress approaches are 
identical. Usually the solution based on the macroscopic strain approach is more simply 
implemented than the macroscopic stress approach. This is because in the macroscopic 
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strain approach, the displacement at the boundary can be simply obtained using Eqs. (5 or 
6). In the macroscopic stress approach, one needs to determine the forces at each 
boundary node and pay special attention to the comer boundary nodes which lie on 
multiple boundary surfaces; proper constraint must also be introduced to assure 
uniqueness of the solution. 

Table 1 lists the expression for C and B for the orthogonal cell (Fig. la). For this 
cell, there is no bending effect on the compliance components for Si i and S 22 - For S44, the 
cell compliance is determined only by the bending. Obviously, for this structure the 
Poisson’s ratio is zero which leads to zero S 12 . Table 2 shows the parameters for a 
triangular cell (Fig. lb). For this structure, all compliances can be written in the form of 
an extension effect but C is dependent on the bending moment I. Table 3 gives the 
parameters for a general honeycomb cell. It shows that the extension and bending 
coefficients C and B depend on the angles between the stmts 0. 

The foam structures shown in Fig. 1 are formed by only one unit cell, given at the 
bottom of these figures. The compliances listed in Table 1-3 also represent the 
compliances for the corresponding foam structures. To show the accuracy of the 
analytical equations, we calculated the elastic properties for the foam structures with 
27*27 cells by solving the system numerically. Figure 1 1 compares the Young’s modulus 
E 11 and shear modulus G 12 for a cellular structure with square cells (L=L z =L y ) at different 
stmt widths t. It shows that the analytical and numerical results are identical. As 
expected, the Young’s modulus is a linear function of t/L and the shear modulus is 
proportional to (t/L) 3 . Figure 12 compares the Young’s modulus Eu and shear modulus 
G /2 for a cellular structure with triangular cells. The two solutions give almost identical 
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results. Both moduli have nearly linear relation with t/L. Figure 13 presents the results for 
a structure with regular hexagonal cells (0=30°). The analytical and numerical solutions 
are slightly different from each other and the discrepancy increases as t/L increases. It has 
been observed that the numerical solution will converge to the analytical solution with 
increase of number of cell in the cellular structures. Therefore we believe this 
discrepancy is introduced by the disturbing of the struts on the boundaries. Further 
comparison for general hexagonal cells is given in Fig. 14 for cells with different 0. For 
t/L=0.2 the solutions are compared with the analytical equations given in [2, 14]. This 
shows that the analytical compliances given in Table 3 have very good agreement with 
the numerical results calculated using 27*27*27 cells for all angles 6 . When the angles 
become negative, the Poisson’s ratio becomes negative also (the so-called reentrant 
structures [19]). 

It is useful to represent the elastic properties in terms of the density ratio. For low 
density foams (t/L« 1) and neglecting the shear correction constant, one can obtain the 
S\ i and 1 S 44 for square foams as 

S n =—> S 44 = — ■ ( 11 ) 

P P 

^ii, Si 2 and 1 S 44 for triangular foams as 

S..=->S 12 ' (12) 

P P P 


S’] 1 , Si 2 and S 44 for regular honeycomb foams as 


s 2 c = _L__L , 2 J_ 

11 2p 3p’’ " 2p 3p !, “~p 3p ! ' 


03 ) 
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Equations (11, 12 and 13) are identical with those shown in refs. 2 and 4. However, the 
compliance list in tables 1-3 can apply to foams with relative higher density ratio. 

5. Foams with 3-D Microstructures 

As for two-dimensional structures, we apply the algorithm described in section 3 
for 3-D orthogonal, rhombic dodecahedron and truncated octohedron structures shown in 
Fig. 2. The unit cell for these three structures are shown at the bottoms of Fig. 2. The 
compliances obtained for these unit cells are listed in Table 4-6. Similarly to the two- 
dimensional compliance matrices, the three-dimensional compliances also can be 
decomposed to axial (C) and bending (B) compliances and the results obtained from the 
macroscopic strain and stress approaches are identical. The Voronoi unit cells for 
rhombic dodecahedron and truncated octahedron are not the smallest repeating cells of 
the corresponding structures. From the structure of each node, one can find the smallest 
repeating cells which are shown in Figure 16 (a) for rhombic dodecahedron and (b) 
truncated octahedron. The effective elastic properties of the foams with truncated 
octahedron cell structures have been studied based on stress analysis using the structural 
beam model [16, 17]. 

Table 4 gives the formulations of C and B for a 3-D orthogonal cell shown of the 
bottom of Fig. 2(a). For this cell there is no bending compliance for S„ (/- 1,2,3), while 
for the shear compliances (S44, S55 and S55) are determined only by bending. Obviously, 
for this structure the Poisson’s ratio is zero which leads to zero S 12 . For the rectangular- 
parallelepiped cellular structures shown at the top of Fig. 2(a), there is only the 
orthogonal cell, therefore the cell compliances shown in Table 4 also represent the 
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compliance for the cellular structures. We compared it with numerical results calculated 
using 8*8*8 cells. The comparison is shown in Fig. 15. 

Table 5 shows the compliances for the smallest periodic cell (Fig. 3(b)) of the 
cellular structures with rhombic dodecahedron cells. The rhombic dodecahedron cellular 
structure is formed by the cell given in Fig. 3(b) with a same-size vacuum. Therefore the 
compliances for this cellular structure will be twice the compliances given in Table 5. 
The comparison with numerical results calculated using 8*8*8 cells are shown in Fig. 17. 

Table 6 shows the compliances for the smallest periodic cell (Fig. 2(c)) of the 
cellular structures with truncated octahedron cells. As shown in Figure 2(c), each node 
has the same structure as the smallest cell shown in Fig. 16(b). However half of the nodes 
whose corresponding structures is rotated by 90 degree from the cell given in Fig. 16(b). 
Therefore the compliances S £ and Poisson’s ratio v,j for this cellular material are given 
by 


=(S C (1,1) + S C (3,3)) / 2 , 5 £ =(5 c (4,4) + 5 c (5,5))/2, 

, 5 ^ 03 ) 

13 S c (3,3) ’ 


(14) 


Comparison with numerical results is given in Fig. 18. 

As in the2D case, one may represent the elastic properties in terms of density ratio for 
low-density foams (t/L« 1) by neglecting the shear correction constant. 


6. Effect of randomness 

All practical foams have a certain degree of randomness. The effects of non- 
periodic microstructure on the elastic properties of 2D and 3D foams are studied by using 


19 



a “drop” nuclei algorithm and Voronoi tessellation to generate foams with different 
regularities [8, 14, 28]. The effect of microstructural heterogeneities on the elastic 
modulus and yield strength has also been investigated by introducing random 
perturbations into the regular honeycomb foams [29]. All these simulations show that 
randomness introduces modulus increases. The aim of this section is to investigate 
quantitatively the relation between elastic properties and packing fraction using the 
algorithm described in section 2, where we showed that the packing fraction is a good 
parameter to describe the regularity of a random structure. At low regularity, the 
geometric properties (size, shape) vastly differ among the cells in the foam. As the 
packing fraction increases, the cells become more uniform and resemble the hexagon in 
2D and the truncated octahedron in 3D. Therefore packing fraction may be used as a 
scaling factor to relate the elastic properties between random and regular foams. 

6. 1 Elasticity of 2D random foams 

Using the microstructure simulation algorithm described in section 2, we obtain 
1 8 sets of independent microstructures for each packing density. Replacing each edge by 
a uniform strut we obtain foam structures with different regularities. Here we use 
constant density ratio for all random structures to avoid the density ratio decrease with 
packing fraction when the cross section area of the strut is a constant. Therefore we vary 
the cross section area to keep the density ratio constant; i. e. the cross section area 
increases when the packing fraction increases. Each edge is modeled by one beam 
element and using the algorithm described in section 3, we obtain the elastic properties of 
the foams. 
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Figures 19(a, b, c) show the Young’s modulus, shear modulus and Poisson ratio 
variation as packing density increases. The density ratio of the foam is 0.1. Each solid 
line represents one independent data set. The average value is shown by the thick solid 
line with open circles. The Young’s modulus and shear modulus decrease as packing 
density increases. The change from the fully random state to the random closed packing 
state is about 30% for both Young’s modus and shear modulus. The Poisson ratio slightly 
increases as packing fraction increases (only about 3% change). The scattering of the data 
is greater at lower packing fractions. 

In Figure 20 we show the modulus as a function of density ratio for different 
packing fraction. Comparing this figure with Fig. 13 for honeycomb structure, one sees 
that the relation between modulus and density ratio is similar. In section 3, we show that 
the compliances of a given structure can be expressed as a summation of bend and axis 
compliances. For random structures, we may also decompose the compliance into 
bending and axis deformation effects 



(15) 


The parameters D % , B E , D G R and B R can be found by fitting the curves given in Figure 
20. These parameters are functions of foam geometrical properties which can be 
categorized by the packing fraction. Figure 20(a) shows the fitting of the data obtained by 
the FEM method (open symbols) using Eq. (15) (solid lines). The parameters obtained 
are shown in Figure 21. The bending parameters B E R and B G R increase as packing fraction 
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increases while the axis deformation parameters D £ and decrease. Their average 
value may be approximated as linear functions of the packing fraction: 

^| = 1.18(<h-i) + l ,||- = -0.341(i,-<i> r ) + l (16) 

V h *>h 

^- = 3.088(^-^ r ) + l , || = -0.277(<|> A -<|> r ) + l (17) 

U h £>h 

where Dj; (=3/2), Bf; (=2/3), D° (=2) and B° (=8/3) are the corresponding parameters for 

the regular honeycomb structure. Therefore, using the packing fraction, we statistically 
correlated the elastic properties for the regular honeycomb to fully random structures. 

6.2 Elasticity of 2D random foams 

As in the 2D case, we simulated 18 sets of independent microstructure for each 
packing density using the algorithm described in Section 2. Each edge is modeled by a 
beam element and using the algorithm described in section 3, we obtain the elastic 
properties of the foams. 

Figures 22(a, b, c) show the Young’s modulus, shear modulus and Poisson ratio 
change as packing density increases. The density ratio of the foam is 0.1. As in the 2D 
case the Young’s modulus and shear modulus decrease as packing density increases and 
the Poisson ratio slightly increases. Because the cell number is smaller, the scattering is 
much greater in the 3D case. 

7. Conclusion 
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In this report analytical and numerical methods to estimate the effective elastic 
properties of regular and random open-cell foams have been developed. The methods are 
based on the principle of minimum energy and on structural beam models. The analytical 
solutions are obtained by symbolic solution of the system. The microstructure of the 
foams is simulated using Voronoi tessellation. By comparing with solid finite element 
solution, the applicability of the solution based on the beam model has been discussed. 
The effect of randomness is also discussed. 
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Figure 1. 2-D microstructures, (b) Rectangle, (b) Triangle, (c) General honeycomb. 


29 









Figure 4. Evolution of 2-D foam microstructure as packing fraction Fc increases. 
Simulated using 1000 random points in a unit square area. 
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Figure 5. Evolution of 3-D foam microstructure as packing fraction Fc increases. 
Simulated using 1000 random points in a unit square area. 
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Figure 6. Statistics of 2D random Voronoi diagram for different regularities, (a) strut 
length distribution, (b) cell perimeter length, (c) cell area distribution, (d) strut angle 
distribution, (e) number of edges per cell. 
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Packing fraction 


Figure 7. 2-D foam density ratio versus random packing fraction for constant strut cross 
section. Thin solid lines correspond to independent random structures at different packing 
fraction. The open circles represent the average, the thick solid lines the linear 
approximation. 
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Figure 8. Statistics of 3D random Voronoi diagram for different regularities, (a) edge 
length distribution, (b) face area distribution, (c) vertex angle distribution, (d) cell face 
area distribution, (e) cell volume length, (f) number of faces per cell, (g) number of edges 
per face. 
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Packing fraction 


Figure 9. 3-D foam density ration versus random packing fraction. Nine sets of random 
structures are calculated at each packing fraction. The thick line represents the average. 
The cross section of the strut is 0.0024 
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Table 1. Cell compliances for general rectangular 2-D 
structures. 
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Figure 11. Comparison of elastic properties calculated numerically using 8*8*8 cells and 
the analytical results listed in Table 1. 
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Table 2. Cell compliances for general triangular 2-D structures. 
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Figure 12. Comparison of elastic properties calculated numerically using 8*8*8 cells and 
the analytical results listed in Table 2. 
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Table 3. Cell compliances for general hexagon 2-D structures. 
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Figure 13. Comparison of elastic properties calculated numerically using 8*8*8 cells and 
the analytical results listed in Table 3. 






Figure 14. Comparison of elastic properties calculated using finite element beam method 
and analytical solutions for general hexagonal structure with different angles 9 . 
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Table 4. Cell compliances for 3-D cubic structures. 
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Figure 15. Comparison of elastic properties calculated numerically using 8*8*8 cells and 
the analytical results listed in Table 4. 
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Table 5. Cell compliances for rhombic dodecahedron structures. 
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Figure 17. Comparison of elastic properties calculated numerically using 8*8*8 cells and 
the analytical results listed in Table 5. 
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Table 6. Cell compliances for truncated octahedron 3D structures. 
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Figure 18. Comparison of elastic properties calculated numerically using 8*8*8 cells and 
the analytical results listed in Table 6. 
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Figure 19. 2-D Elastic properties versus random packing fraction. 18 sets of random 
structures are calculated at each packing fractures. The thick line represents the average. 
The density ratio of the foam is 0.1. 
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Figure 20. 2-D Elastic properties versus relative density for different packing fractions. 
18 sets of random structures are calculated at each packing fraction. 
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Figure 21 . 2-D elastic properties versus random packing fraction. Nine sets of random 
structures are calculated at each packing fraction. The thick line represents the average. 
The density ratio of the foam is 0.1. 
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Packing fraction 


Figure 22. 3-D elastic properties versus random packing fraction. Nine sets of random 
structures are calculated at each packing fraction. The thick line represents the average. 
The density ratio of the foam is 0. 1 . 
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